Satellite Image Data


Analysis using numpy


Data Source: Satellite Image from WIFIRE Project

WIFIRE is an integrated system for wildfire analysis, with specific regard to changing urban dynamics and climate. The system integrates networked observations such as heterogeneous satellite data and real-time remote sensor data, with computational techniques in signal processing, visualization, modeling, and data assimilation to provide a scalable method to monitor such phenomena as weather patterns that can help predict a wildfire's rate of spread.

In this example, we will analyze a sample satellite image dataset from WIFIRE using the numpy Library.

Loading the libraries we need: numpy, scipy, matplotlib

In [2]:
%matplotlib inline
import numpy as np
from scipy import misc
import matplotlib.pyplot as plt

Creating a numpy array from an image file:


Lets choose a WIFIRE satellite image file as an ndarray and display its type.

Download the image from here

In [3]:
#from skimage import data

photo_data = misc.imread('./wifire/sd-3layers.jpg')

type(photo_data)
Out[3]:
numpy.ndarray

Let's see what is in this image.

In [4]:
plt.figure(figsize=(15,15))
plt.imshow(photo_data)
Out[4]:
<matplotlib.image.AxesImage at 0x29896975240>
In [5]:
photo_data.shape

#print(photo_data)
Out[5]:
(3725, 4797, 3)

The shape of the ndarray show that it is a three layered matrix. The first two numbers here are length and width, and the third number (i.e. 3) is for three layers: Red, Green and Blue.

RGB Color Mapping in the Photo:


  • RED pixel indicates Altitude

  • BLUE pixel indicates Aspect

  • GREEN pixel indicates Slope

    </ul>
    The higher values denote higher altitude, aspect and slope.

In [10]:
photo_data.size
Out[10]:
53606475
In [11]:
photo_data.min(), photo_data.max()
Out[11]:
(0, 255)
In [17]:
photo_data.mean()
Out[17]:
75.829935450894695


Pixel on the 150th Row and 250th Column

In [19]:
photo_data[150, 250]
Out[19]:
array([ 17,  35, 255], dtype=uint8)
In [18]:
photo_data[150, 250, 1]
Out[18]:
35


Set a Pixel to All Zeros


We can set all three layer in a pixel as once by assigning zero globally to that (row,column) pairing. However, setting one pixel to zero is not noticeable.

In [20]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
photo_data[150, 250] = 0
plt.figure(figsize=(10,10))
plt.imshow(photo_data)
Out[20]:
<matplotlib.image.AxesImage at 0x29896b782e8>


Changing colors in a Range


We can also use a range to change the pixel values. As an example, let's set the green layer for rows 200 t0 800 to full intensity.

In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')

photo_data[200:800, : ,1] = 255
plt.figure(figsize=(10,10))
plt.imshow(photo_data)
In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')

photo_data[200:800, :] = 255
plt.figure(figsize=(10,10))
plt.imshow(photo_data)
In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')

photo_data[200:800, :] = 0
plt.figure(figsize=(10,10))
plt.imshow(photo_data)


Pick all Pixels with Low Values

In [22]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
print("Shape of photo_data:", photo_data.shape)
low_value_filter = photo_data < 200
print("Shape of low_value_filter:", low_value_filter.shape)

Shape of photo_data: (3725, 4797, 3)
Shape of low_value_filter: (3725, 4797, 3)

Filtering Out Low Values


Whenever the low_value_filter is True, set value to 0.

In [23]:
#import random
plt.figure(figsize=(10,10))
plt.imshow(photo_data)
photo_data[low_value_filter] = 0
plt.figure(figsize=(10,10))
plt.imshow(photo_data)
Out[23]:
<matplotlib.image.AxesImage at 0x29896f2b5c0>

More Row and Column Operations


You can design complex patters by making cols a function of rows or vice-versa. Here we try a linear relationship between rows and columns.

In [27]:
rows_range = np.arange(len(photo_data))
cols_range = rows_range
print(type(rows_range))

<class 'numpy.ndarray'>
In [28]:
photo_data[rows_range, cols_range] = 255
In [29]:
plt.figure(figsize=(15,15))
plt.imshow(photo_data)
Out[29]:
<matplotlib.image.AxesImage at 0x29898256a90>


Masking Images


Now let us try something even cooler...a mask that is in shape of a circular disc.

In [51]:
total_rows, total_cols, total_layers = photo_data.shape
#print("photo_data = ", photo_data.shape)

X, Y = np.ogrid[:total_rows, :total_cols]
print("X = ", X.shape, " and Y = ", Y.shape)
X

X =  (3725, 1)  and Y =  (1, 4797)
Out[51]:
array([[   0],
       [   1],
       [   2],
       ..., 
       [3722],
       [3723],
       [3724]])
In [49]:
center_row, center_col = total_rows / 2, total_cols / 2
#print("center_row = ", center_row, "AND center_col = ", center_col)
#print(X - center_row)
#print(Y - center_col)
dist_from_center = (X - center_row)**2 + (Y - center_col)**2
print(dist_from_center.shape)
print(dist_from_center[:10,:10])
radius = (total_rows / 2)**2
#print("Radius = ", radius)
circular_mask = (dist_from_center > radius)
#print(circular_mask.shape)
#print(circular_mask)
#print(circular_mask[1500:1700,2000:2200])

(3725, 4797)
[[ 9221708.5  9216912.5  9212118.5  9207326.5  9202536.5  9197748.5
   9192962.5  9188178.5  9183396.5  9178616.5]
 [ 9217984.5  9213188.5  9208394.5  9203602.5  9198812.5  9194024.5
   9189238.5  9184454.5  9179672.5  9174892.5]
 [ 9214262.5  9209466.5  9204672.5  9199880.5  9195090.5  9190302.5
   9185516.5  9180732.5  9175950.5  9171170.5]
 [ 9210542.5  9205746.5  9200952.5  9196160.5  9191370.5  9186582.5
   9181796.5  9177012.5  9172230.5  9167450.5]
 [ 9206824.5  9202028.5  9197234.5  9192442.5  9187652.5  9182864.5
   9178078.5  9173294.5  9168512.5  9163732.5]
 [ 9203108.5  9198312.5  9193518.5  9188726.5  9183936.5  9179148.5
   9174362.5  9169578.5  9164796.5  9160016.5]
 [ 9199394.5  9194598.5  9189804.5  9185012.5  9180222.5  9175434.5
   9170648.5  9165864.5  9161082.5  9156302.5]
 [ 9195682.5  9190886.5  9186092.5  9181300.5  9176510.5  9171722.5
   9166936.5  9162152.5  9157370.5  9152590.5]
 [ 9191972.5  9187176.5  9182382.5  9177590.5  9172800.5  9168012.5
   9163226.5  9158442.5  9153660.5  9148880.5]
 [ 9188264.5  9183468.5  9178674.5  9173882.5  9169092.5  9164304.5
   9159518.5  9154734.5  9149952.5  9145172.5]]
In [41]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
photo_data[circular_mask] = 0
plt.figure(figsize=(15,15))
plt.imshow(photo_data)
Out[41]:
<matplotlib.image.AxesImage at 0x298982fb668>

Further Masking


You can further improve the mask, for example just get upper half disc.

In [ ]:
X, Y = np.ogrid[:total_rows, :total_cols]
half_upper = X < center_row # this line generates a mask for all rows above the center

half_upper_mask = np.logical_and(half_upper, circular_mask)
In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
photo_data[half_upper_mask] = 255
#photo_data[half_upper_mask] = random.randint(200,255)
plt.figure(figsize=(15,15))
plt.imshow(photo_data)


Further Processing of our Satellite Imagery

Processing of RED Pixels

Remember that red pixels tell us about the height. Let us try to highlight all the high altitude areas. We will do this by detecting high intensity RED Pixels and muting down other areas.

In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
red_mask   = photo_data[:, : ,0] < 150

photo_data[red_mask] = 0
plt.figure(figsize=(15,15))
plt.imshow(photo_data)


Detecting Highl-GREEN Pixels

In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
green_mask = photo_data[:, : ,1] < 150

photo_data[green_mask] = 0
plt.figure(figsize=(15,15))
plt.imshow(photo_data)


Detecting Highly-BLUE Pixels

In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')
blue_mask  = photo_data[:, : ,2] < 150

photo_data[blue_mask] = 0
plt.figure(figsize=(15,15))
plt.imshow(photo_data)


Composite mask that takes thresholds on all three layers: RED, GREEN, BLUE

In [ ]:
photo_data = misc.imread('./wifire/sd-3layers.jpg')

red_mask   = photo_data[:, : ,0] < 150
green_mask = photo_data[:, : ,1] > 100
blue_mask  = photo_data[:, : ,2] < 100

final_mask = np.logical_and(red_mask, green_mask, blue_mask)
photo_data[final_mask] = 0
plt.figure(figsize=(15,15))
plt.imshow(photo_data)
In [ ]: